Novel Method of Measuring Corneal Viscoelasticity Using the Corvis ST Tonometer

Background: The process of rapid propagation of the corneal deformation in air puff tonometer depends not only on intraocular pressure, but also on the biomechanical properties of the cornea and anterior eye. One of the biomechanical properties of the cornea is viscoelasticity, which is the most visible in its high-speed deformations. It seems reasonable to link the corneal viscoelasticity parameter to two moments of the highest speed of corneal deformations, when the cornea buckles. The aim of this work is to present a method of determining the time and place of occurrence of corneal buckling, examine spatial and temporal dependencies between two corneal applanations and bucklings in the Corvis ST tonometer, and correlate these dependencies with corneal viscoelastic properties. Methods: Images of the horizontal cross section of the Corvis ST deformed cornea from the air puff tonometer Corvis ST were used. 14 volunteers participated in the study, each of them had one eye measured eight times. Mutual changes in the profile slopes of the deformed corneas were numerically determined. They describe pure corneal deformation, eliminating the influence of rotation, and displacement of the entire eyeball. For each point in the central area of the corneal profile, the maximum velocities of mutual slope changes accompanying the applanations were estimated. The times of their occurrence were adopted as buckling times. Results: The propagation of buckling along the corneal profile is presented, as well as the repeatability and mutual correlations between the buckling parameters and intraocular pressure. Based on the relationship between them, a new parameter describing corneal hysteresis: Corvis Viscoelasticity (CVE) is introduced. It is characterized by high repeatability: ICC = 0.82 (0.69–0.93 CI) and low and insignificant correlation with intraocular pressure: r = 0.25 (p-value = 0.38). Conclusion: The results show for the first time how to measure the corneal buckling and viscoelastic effects with Corvis ST. CVE is a new proposed biomechanical parameter related to the viscoelastic properties of the cornea, which has high repeatability for the examined subject. The distribution of its values is planned to be tested on different groups of patients in order to investigate its clinical applicability.


Introduction
Variations in the mutual inclinations of the tangents to the corneal profile are one way to describe corneal deformation. During air puff tonometric measurements, these mutual inclinations are subjected to change, conveying information about pure deformation of the cornea, being insensitive to shifts or rotation of the entire eyeball caused by the air blast. Corneal deformation depends both on the resultant forces that act on it and its biomechanical properties.
The Corvis ST tonometer (OCULUS Optikgeräte GmbH, Wetzlar, Germany) records a sequence of images during each measurement, showing corneal deformations in the horizontal cross-section induced by an air blast. These images can be used to compute the described above mutual slopes of the various corneal regions.
The device provides information on intraocular pressure IOP, thickness distribution of the horizontal corneal profile, and a number of dynamic corneal response (DCR) parameters that describe the behavior of the cornea during measurement. Taking into account the biomechanical properties of the cornea can be useful as an additional diagnostic tool. For example, the parameters for the dynamic corneal response may indicate the early stage of keratoconus even when tomography and topography appear normal [1]. Furthermore, in normal-tension glaucoma, corneas tend to be softer and more deformable compared to the control group, which is reflected in their dynamic response to air puff [2].
The repeatability of individual parameters and their mutual relationships are important aspects of data analysis. Particular attention is paid to the relationship between the DCR parameters and the intraocular pressure. This is because the IOP is measured indirectly through the cornea. Its biomechanical properties (thickness, stiffness, and viscoelasticity) can significantly distort the reading. It also works the other way around: due to stresses induced in the viscoelastic cornea by varying intraocular pressure, its biomechanical parameters may change depending on the current IOP value [3]. This issue is extremely complex and new models are needed to address it.
The biomechanically corrected intraocular pressure bIOP is a characteristic parameter determined by this tonometer. The equation for the bIOP value was determined using numerical simulations with the finite element method, taking into account the influence of corneal stiffness, corneal thickness, curvature, and biomechanical properties [4]. It has been shown experimentally that it is in fact less dependent on corneal thickness and age than uncorrected IOP, measured with the same tonometer [5,6]. Recently, a material stiffness parameter (Stress-Strain Index SSI) has also been proposed, which manifests a high correlation with age and a lack of significant correlation with IOP and corneal thickness. Based on the image sequence from the Corvis ST tonometer, SSI allows one to estimate the tangent modulus of the cornea under any load or stress [7].
The parameters proposed by the manufacturer are just some of many possible ways of describing, using, or interpreting corneal deflection. However, the complex phenomena that occur during full corneal deflection are not yet fully understood. There is a need to learn more about them, as they can tell more about the properties of the cornea itself, as well as about the entire eye globe or the structures surrounding it.
The corneal stroma consists of about 78% water [8]. Fast deformations of structures containing liquids are accompanied by viscous resistance inside this structure. Viscous resistance depends directly on the velocity of deformation in such a structure [9]. The viscous resistance of cornea can be expected to depend on the velocity of the corneal deformation. During air puff measurements, the fastest deformations inside the cornea appear approximately when the central part of the cornea passes both applanations. This suggests that the corneal viscoelasticity directly affects the way the cornea passes both applanations.
Viscoelastic properties have not yet been described or examined with the Corvis ST tonometer. However, they are determined in the Ocular Response Analyzer (ORA), using corneal hysteresis CH, which is linearly dependent on the difference between the air pressures P 1 and P 2 from the jet, during the first and the second applanations [10]. CH is used to diagnose and manage suspected glaucoma. It is also useful for monitoring changes, especially in the case of normal tension glaucoma [11]. Some of the recently released articles constate that the viscous properties of the cornea cannot be defined from the air puff applanation. One of these articles suggests that the viscous properties of the cornea did not contribute significantly to the measured displacements [12], but another states that any fast contactless test can provide information on the viscous properties of the cornea [13]. Moreover, Francis et al. directly argue that corneal viscous properties cannot be determined from air puff tonometry [14]. However, this article was commented on by other authors, who stated that hysteresis-an indicator of viscoelasticity-can be determined using both Corvis ST and ORA [15].
Usually, when considering the material viscoelasticity, it is assumed that the material is linearly elastic. A more general approach describes the hyperelastic-viscoelastic material.
Such materials are often known as nonlinear viscoelastic materials. It is mentioned in the literature that some biological tissues indicate these hyper-viscoelastic properties. Zhang et al. consider both viscoelasticity and hyperelasticity to establish a more accurate model of the human ear [16]. This approach was also proposed for the material properties of the cornea. Jannesari et al. state that the contribution of material viscosity in the dynamic tonometry test is trivial and can be ignored, but hyper-viscoelastic is a more accurate definition of the behavior of the whole eye globe [17]. The viscoelastic-hyperelastic model of the cornea is examined in papers [18,19]. However, examinations of both models were carried out on ex vivo porcine and human corneas, during relatively slow, monotonic loading, and deformations. Corneal deformations during real-time air puff tonometry are significantly faster and distinctly nonmonotonic.
The process of very fast corneal applanation is very complex and is not well recognized from a mechanical point of view. Theoretical and experimental analysis of the complex structure, quasi-spherical, viscoelastic shell buckling would be essential and very interesting for modeling and investigating corneal buckling during its air puff applanation. Some papers concern buckling of the elastic spherical shell [20]. However, to the best of our knowledge, there are no such papers concerning the viscoelastic shell.
In the paper [21], both corneal applanations in ORA are examined and the corneal buckling process is described and discussed. It draws attention that the times referring to both bucklings t A and t B are very close to both applanation times t A1 and t A2 , but they are not precisely the same.
Numerical analysis of the dynamics of the corneal profile inclinations observed in image sequences from Corvis ST may enable the identification of two-time moments, when the two mutual inclinations of the corneal profile in the area of the apex are extremely fast. If these time moments are very close to times of the first and the second applanations, they can be assumed to be the times of both corneal bucklings. Since buckling of the central cornea due to air blast is associated with a sudden and rapid redistribution of stresses inside the cornea, it is likely that it can also be observed in Corvis ST measurements during both applanations.
The aim of this work is to present the method for determining the time and place of occurrence of corneal buckling, to examine spatial and temporal dependencies between the two corneal applanations and the bucklings in the Corvis ST tonometer, and to correlate these dependencies with corneal viscoelastic properties. The proposed parameters may become an alternative to corneal hysteresis measured by ORA, achievable in the measurements with Corvis ST.

Materials and Methods
These retrospective data were collected during a previous study, conducted at an outpatient unit of the Department of Ophthalmology of the Medical University Hospital in Wroclaw, Poland, with the prior consent of the bioethics committee. We did not perform the sample size calculation, as this is the pilot study that we intend to use in the future. Based on this perspective, it was deemed important to assess the possible best-choice parameters on the smaller group, as the analytical approach is a time-and resourceconsuming procedure that requires several iterations of the process. Fourteen healthy volunteers participated in this study (the average age of the participants was 28 years with a standard deviation of ±11). The tested group included healthy Caucasian men and women, without diagnosed eye diseases, symptoms, or previous eye surgery. They were examined using the noncontact Corvis ST tonometer (OCULUS Optikgeräte GmbH). During each measurement, the tonometer was set in the appropriate position just in front of the subject's eye, according to the instructions displayed on the device screen. A blast of air was released from the device nozzle to deflect the cornea.
Each person had performed eight measurements on one selected eye, on the same day, one after the other, several minutes apart. Subjects were encouraged to relax and move their heads between subsequent measurements. The whole procedure for a single eye lasted about 30 min. These multiple measurements were taken on seven right eyes and seven left eyes. The allocation of people to the left or right eye measurement group was random (ocular dominance was not taken into account). All subjects were informed about related procedures and consented to the examination. The study was carried out according to the Declaration of Helsinki and obtained ethical approval.
During each measurement, the device allowed the registration of a sequence of 140 images (200 × 576 pixels) displaying a horizontal cross-section of a deforming cornea. Each sequence lasted about 32 ms. The size of a pixel was assumed to be 0.016 × 0.016 mm [22]. A total of 112 image sequences were exported using the Corvis ST software (version 1.3r1727) and processed using Matlab.
First, horizontal profiles of the examined corneas were detected, as described in detail in [23]. As a result, for each measurement, two-dimensional matrices M(i, k) were obtained (number of images i from 1 to 140, which corresponds to the period from 0 to 32 ms; image columns k from 1 to 576). These matrices were smoothed using a 5 × 7 window size averaging filter to eliminate possible disturbances. This smoothing procedure created a matrix of smoothed data Ms(i, k), describing the height of the corneal profile in image number i and image column k.
Local slopes α of the corneal profiles were calculated as follows: where d dk is a five-point numerical derivative of the coordinate (column) k. Slope changes, which means slope variation relative to the initial state in the first frame, were calculated according to the following formula: To remove the influence of possible eye rotation on the analyzed values of the slope along the corneal profile, the change in the mutual slope at two points located symmetrically to the apex: Profile Slope Variations (PSV), were determined as follows Since the largest changes due to air blast are observed in the central part of the cornea, our further analysis focused on this area. Two points located approximately 0.768 mm from the center on the nasal and temporal side were selected (corneal profile points falling on the 240th and 337th columns of images-see Figure 1). This choice has been made due to the high reproducibility and correlations of the proposed parameters referred to these points. Such a distance from the apex is large enough for clear observation of the mutual slope changes (low risk of numerical errors caused by too small mutual slope changes), yet much less than the minimum applanation lengths or peak distance along the corneal profile indicated by the tonometer software. Figure 2 shows an example of a PSV function during measurement. In the time moments close to the applanation, PSV dynamically changes its value, while in the remaining moments of the measurement, it fluctuates around certain values.
After a tenfold densification of the PSV function points, with the smoothing spline method, the PSV time derivative was calculated. Two moments B1T and B2T close to the applanation times, in which the velocity (derivative) of the PSV reaches two extreme values DSV1 and DSV2 were determined. Since times B1T and B2T are characterized by the fastest possible corneal deformations (the fastest changes in the local corneal slope), one can assume that they correspond to the times when the complex corneal structure loses its stability under extremely fast loading and starts its buckling.
The times of applanations A and buckling B are very close to each other, but they do not coincide (as shown in Figure 3). The corresponding differences in both times are marked as D1T and D2T: . Clin. Med. 2021, 10, x FOR PEER REVIEW   After a tenfold densification of the PSV function points, with the sm method, the PSV time derivative was calculated. Two moments B1T and B applanation times, in which the velocity (derivative) of the PSV reaches tw   After a tenfold densification of the PSV function points, with the sm method, the PSV time derivative was calculated. Two moments B1T and B applanation times, in which the velocity (derivative) of the PSV reaches tw ues DSV1 and DSV2 were determined. Since times B1T and B2T are chara fastest possible corneal deformations (the fastest changes in the local corn can assume that they correspond to the times when the complex corneal its stability under extremely fast loading and starts its buckling.
The times of applanations A and buckling B are very close to each oth A further step was the correlation analysis between the para and with the parameters obtained from the Corvis ST so noncorrected intraocular pressure, bIOP-biomechanically pressure, CCT-central corneal thickness, A1T and A2T-times applanation, A1L and A2L-lengths of the first and the second a

Statistical Analysis
The ranges and repeatability of the basic parameters obtai the corneal slope parameters were examined. For this purpose, were determined for all eyes measured: mean values, standard d correlation coefficient ICC [24] coefficient of repeatability CR [25 CV. In addition, the incidence of significant relationships betw was also estimated. Based on the mean values for individual coefficients r for pairs of different variables were calculated [26]. level was 0.05. Figure 5 present the graphs of dependencies b well as between D2T and B2T from one measurement as a functi the corneal center. The courses of the first buckling B1T and the se show that they appear at different times depending on the di center. Analysis of these dependencies was carried out before mo In case of the distance closer to the corneal center, the slope var velocity DSV1 and DSV2 are ignorable. On the other hand, when the analyzed area is no longer related to the deflected corneal ce Additionally, the parameters analyzed for each measurement were ( Figure 3):

Figure 4 and
• T mp -time of the maximum air pressure from the nozzle. This parameter was determined on the basis of air pressure diagrams exported from the tonometer software.
As the air pressure curves were very similar to each other, T mp differed only slightly between measurements and was in the range of 16.23-16.72 ms; • dT1 and dT2-time periods from B1T to T mp and from T mp to B2T; • dAT1 and dAT2-time differences T mp -A1T and A2T-T mp .
Some algebraic operations between buckling parameters represent the asymmetry in the dynamical response of the cornea during its loading and unloading, which is related to the hysteresis of the corneal response. This paper presents the results of the analysis for the following parameters: . A further step was the correlation analysis between the parameters introduced above and with the parameters obtained from the Corvis ST software, such as IOP-noncorrected intraocular pressure, bIOP-biomechanically corrected intraocular pressure, CCT-central corneal thickness, A1T and A2T-times of the first and second applanation, A1L and A2L-lengths of the first and the second applanation.

Statistical Analysis
The ranges and repeatability of the basic parameters obtained from Corvis ST and the corneal slope parameters were examined. For this purpose, the following quantifiers were determined for all eyes measured: mean values, standard deviations, and intraclass correlation coefficient ICC [24] coefficient of repeatability CR [25], coefficient of variation CV. In addition, the incidence of significant relationships between selected parameters was also estimated. Based on the mean values for individual eyes, the 's correlation coefficients r for pairs of different variables were calculated [26]. The adopted significance level was 0.05.

Figures 4 and 5 present the graphs of dependencies between D1T
and B1T, as well as between D2T and B2T from one measurement as a function of the distance d from the corneal center. The courses of the first buckling B1T and the second B2T in both figures show that they appear at different times depending on the distance from the corneal center. Analysis of these dependencies was carried out before more extensive calculations. In case of the distance closer to the corneal center, the slope variations and their angular velocity DSV1 and DSV2 are ignorable. On the other hand, when the distance d is too high, the analyzed area is no longer related to the deflected corneal center. Therefore, the time difference D1T becomes very small and then negative. In this situation, obtained results characterize a low repeatability of the analyzed parameters. It should be noted that it is possible to find a distance d that is large enough for B1T to occur later than A1T (d comparable to the length of the first application or greater). Furthermore, B2T always occurs before A2T, regardless of the distance of the analyzed point of the cornea from the apex.
lin. Med. 2021, 10, x FOR PEER REVIEW difference D1T becomes very small and then negative. In this situati characterize a low repeatability of the analyzed parameters. It shoul possible to find a distance d that is large enough for B1T to occur lat parable to the length of the first application or greater). Furthermore before A2T, regardless of the distance of the analyzed point of the cor It is interesting to note that the first buckling starts before the fir corneal center and expands in time radially from the center. In the buckling (Figure 5), it falls in the center of the cornea also before the The dispersal of the first buckling from the corneal center before the slower than the descent of the second buckling to the center of the co ond applanation, while the average time delay between the first bu applanation D1T is smaller than the similar delay D2T during the sec  It is interesting to note that the first buckling starts before the first applanation in the corneal center and expands in time radially from the center. In the case of the second buckling ( Figure 5), it falls in the center of the cornea also before the second applanation. The dispersal of the first buckling from the corneal center before the first applanation is slower than the descent of the second buckling to the center of the cornea before the second applanation, while the average time delay between the first buckling and the first applanation D1T is smaller than the similar delay D2T during the second applanation.
The corneal applanations and the fastest deformations of the corneal profile in the distance applied in the paper (48 px from the apex) appear clearly at different times for both applanations. Interestingly, in each of the analyzed measurements, buckling always occurred right before the corresponding applanation. Detailed descriptive statistics of the parameters considered can be seen in Table 1.  The corneal applanations and the fastest deformations of the co distance applied in the paper (48 px from the apex) appear clearly a both applanations. Interestingly, in each of the analyzed measuremen Analysis of correlations between considered parameters showed interesting and important properties, which enables a new approach to Corvis ST measurements. We decided to split the analyzed correlation coefficients into two characteristic groups of parameters.
The first group contains parameters highly correlated with IOP. They are mainly defined by both applanation times, buckling times, but also by difference D1T − D2T as well as by both extreme speeds of the slope variations DSV1 and DSV2. The mutual correlations between them can be seen in Table 2.
The second group consists of the time B2T and the relations between D1T, D2T, dT1, and dT2. The main property of these parameters is their high mutual correlation related to the viscoelastic properties of the cornea and low correlations with IOP. The mutual correlations between these parameters can be seen in Table 3. An exemplary correlation plot is presented in Figure 6.
Due to the small values of D1T and D2T, the values of some parameter pairs presented in Tables 2 and 3 are very similar to each other, such as A1T and B1T, A2T and B2T, dT1/dT2 and dAT1/dAT2. However, one can see that the time parameters related to both applanations are higher correlated with IOP, while the time parameters referred to both bucklings have better correlations with viscoelastic properties. A good example is the correlation coefficient between A1T and IOP (r = 0.99), while a similar correlation between B1T and IOP is smaller (r = 0.79). On the other hand, the correlation between dT1/dT2 and B2T is higher (r = −0.94) than the correlation between dAT1/dAT2 and A2T (r = −0.78). There are more such observations among the respective data, which show that it is particularly important to take into account the appropriate parameters among similar ones during the calculation of IOP and bIOP or the calculation of parameters referring to viscoelastic properties. Table 1. The mean values of the considered parameters, their standard deviations (SD), ranges, as well as measures of their repeatability in the form of coefficient of repeatability (CR), coefficient of variation (CV), and intraclass correlation coefficient (ICC) with 95% confidence interval (CI); u-the unit of a given parameter. Analysis of the mutual correlations presented in Table 3 enables the proposal of a new parameter, which characterizes the viscoelastic properties of the cornea. We propose to use the parameter D1T dT1 + D2T dT2 , describing asymmetry of corneal behavior during its loading and unloading, using the characteristic relation between buckling and applanation moments. The selected parameter has high repeatability (ICC = 0.82) and its correlations with other parameters from the table are also high. This parameter is also dimensionless, which can be its advantage. To avoid its complex notation, we propose to call this parameter CVE (Corvis ViscoElasticity). The ranges of the CVE values for individual eyes are shown in Figure 7.

Parameters
To analyze the property of two applanation lengths A1L and A2L in the Corvis ST tonometer, their repeatabilities and correlations with other parameters considered were calculated and examined. Repeatability coefficients for both applanation lengths were very low. Several of the parameters analyzed showed a correlation coefficient with A1L greater than 0.7 ( Table 4). The highest correlation appears for the ratio D1T D2T and amounts 0.85. There was a problem to find a high correlation of the second applanation length A2L with any of the parameters considered.  Due to the small values of D1T and D2T, the values of some p sented in Tables 2 and 3 are very similar to each other, such as A1T an dT1/dT2 and dAT1/dAT2. However, one can see that the time param applanations are higher correlated with IOP, while the time parame bucklings have better correlations with viscoelastic properties. A good relation coefficient between A1T and IOP (r = 0.99), while a similar B1T and IOP is smaller (r = 0.79). On the other hand, the correlation b B2T is higher (r = −0.94) than the correlation between dAT1/dAT2 and A are more such observations among the respective data, which show t important to take into account the appropriate parameters among sim calculation of IOP and bIOP or the calculation of parameters referring erties.
Analysis of the mutual correlations presented in Table 3 enabl new parameter, which characterizes the viscoelastic properties of the to use the parameter 1 1 + 2 2 , describing asymmetry of corneal beha ing and unloading, using the characteristic relation between buckli moments. The selected parameter has high repeatability (ICC = 0.82) with other parameters from the table are also high. This parameter is which can be its advantage. To avoid its complex notation, we propos eter CVE (Corvis ViscoElasticity). The ranges of the CVE values for shown in Figure 7.
To analyze the property of two applanation lengths A1L and A

Discussion
The results presented in this paper enable the description and analysis of corneal buckling in measurements using the Corvis ST tonometer. The main areas of interest were processes related to the first and the second applanations, the values and times of extreme angular speeds of the corneal slope variations, and the asymmetry of both applanations in relation to the maximum air pressure from the jet. Extreme speeds of corneal slope variations appear directly before the first and the second applanations, not during these applanations, as expected. Such extreme speeds of variation of the local corneal slope indicate its local buckling. As was presented in both the paper [21] and in this study, buckling effects did not appear directly during both corneal applanations but very close to these flattenings of the cornea. Thus, it is interesting and important to elucidate the reason for this effect and the dependence of the time difference between these processes and the corneal properties.
From a hydrodynamic point of view, a rapid deformation of structures consisting mainly of liquid causes a viscosity effect, which is proportional to the velocity of the deformation or to its power (nonlinear effect). Since the cornea consists mainly of water, one has to expect the prevalence of viscoelastic effect (viscoelastic resistance) during fast corneal deformations, which appear mainly during the time periods of both applanations. The magnitude of this effect depends on the velocity of deformation of the structure but also on the properties of the structure, described by some characteristic constants. This approach suggests correlating the viscoelastic properties of the cornea (characterized by a constant) with the time delay between the corneal applanation and corneal buckling. The magnitude of this delay should be a function of the corneal structure resistance caused by its viscoelasticity. Numerical analysis of the results obtained showed that there are quite high correlations between the selected parameters related to the mentioned time delays (D1T and D2T) and very low correlations with IOP.
An interesting correlation appears between both buckling times B1T and B2T vs. D1T. No correlation between D1T and B1T means that the time of the first buckling B1T (and the time of the first applanation A1T) does not have any influence on their difference D1T. Both A1T and B1T times depend mainly on the IOP behind the cornea. However, unexpectedly high correlations were obtained between the time difference D1T during the first applanation and the second buckling time B2T as well as the second applanation A2T. This indicates that processes between the first buckling and the first applanation determine the times when the second buckling and the second applanation appear. This effect can be explained by the fact that the two intervals D1T and D2T characterize the viscoelastic properties of the cornea, due to rapid processes within these intervals. Such properties also determine time intervals between applanations (r = −0.71 for D1T vs. A2T − A1T) and bucklings (r = −0.70 for D1T vs. B2T − B1T). The sum of times D1T and D2T underlines the viscoelastic properties of the cornea, while their difference reduces these properties being more related to IOP. Here, some far-reaching analogies can be seen with Ocular Response Analyzer (ORA) measurements, where the sum of air pressures during corneal applanations P1 and P2 is proportional to the IOP value, and the difference of both pressures is proportional to Corneal Hysteresis CH [27,28].
Furthermore, it appears that the parameters that describe the asymmetry of both bucklings dT1 dT2 and applanations dAT1 dAT2 in relation to the time of maximum air pressure T mp are somehow correlated with the parameters associated to the corneal viscoelastic properties.
A similar effect appears in the ORA device. According to our earlier examinations, the correlation coefficient of corneal hysteresis CH with the similar asymmetry dAT1 dAT2 of both applanations falls within the range 0.75-0.85. However, in ORA measurements, dAT1 is smaller than dAT2, opposite to Corvis ST measurements.
Some authors affirm that the viscous properties of the cornea do not contribute significantly to the measured displacements in air puff tonometers or even that they cannot be determined [12,14]. They refer to in vitro measurements on a slowly and monotonically deformed cornea. The average viscous properties of the cornea under such slow deformations probably demonstrate a very small influence on the results obtained. However, in the case of rapid variability in deformation during an extremely short period, when buckling and applanation of the cornea appear within a fraction of milliseconds, the viscous effect can be clearly observed and measured. This effect is most likely nonlinear and is clearly observed only during a quick variation of corneal deformation.
The applanation of the cornea means the flattening of the first corneal surface. Taking into account the increase in corneal thickness outside its center, the flat anterior surface of the cornea means that the borders of the deeper-located layers of the corneal structure are still convex. Due to the very complex corneal structure, it is difficult to describe how the mentioned borders pass their flattening, and how an extremely fast buckling propagates inside a complex, multilayered, anisotropic structure. Moreover, corneal applanation cannot indicate an ideally flat area of the central applanation. There are always rapid vibrations or small irregularities on this surface. Therefore, according to the article [21] the question of 'quality of both corneal applanations' is still open. It cannot be excluded that various corneas are characterized by various 'qualities of applanation' or areas of applanation. The presented correlation coefficients with the first applanation length A1L suggest that this parameter can be somehow related to other examined parameters.
One of the limitations of the study is that it was performed on a small sample of healthy eyes, so it does not provide much information about expected values of population parameters. We are aware that the study sample size should be larger to increase the validity of the new proposed method. Now, with the direction of further research already established, we plan to test the method in larger groups of subjects, both healthy and suffering from some specific pathologies. On the other hand, the average values of the parameters read from the tonometer software and their repeatability are comparable to the values presented in studies using larger study groups [29,30].
Another limitation was that the Corvis ST tonometer cannot change the characteristics of the air pulse-its length and pressure. These parameters are constant and repeatable in subsequent measurements, and therefore it is not possible to examine their influence on the considered parameters. In addition, it would be good to model the buckling in the cornea and adjust the mechanical parameters of the model to reflect its actual behavior. However, until an exact model is not available, this step cannot be accomplished.
The other possible drawback of the study is the fact that we based only on the commercially available devices. However, scientific research is carried out on OCT devices to record the reaction of eye structures to a blast of air [31]. Possibly, the parameters related to corneal buckling could also be determined with this newly developed device.
Therefore, in order to unify the measurement results obtained on each of these devices, their different air impulse characteristics as well as different methods of recording the dynamic reaction of the cornea (corneal reflection light detector, Scheimpflug images or OCT) should be taken into account.
During the numerical analysis of the repeatability of the parameters and their mutual correlations, the influence of the Central Corneal Thickness CCT on the values obtained was also examined. The parameters were multiplied or divided by CCT and then analyzed. In some cases, such an operation increases the parameter repeatability, and in other cases increases their mutual correlation. However, this influence of CCT on the results received was not very clear and unique. More detailed analysis of this effect could be carried out in the next examinations.
The use of Corvis ST to measure dynamic changes of the cornea allows for a noninvasive study of the biomechanical properties of the eye. We believe that the new parameters introduced, which show a high repeatability and very low correlation with intraocular pressure, are highly correlated with the viscoelastic properties of the cornea during applanation. The proposed CVE parameter may find application in the diagnosis of the early stage of corneal keratoconus, the variation of corneal properties due to glaucoma, the effect of corneal crosslinking, or the influence of corneal viscoelastic properties on different refractive surgery procedures.
CVE describes the viscoelastic properties of the cornea, in a way similar to the CH hysteresis in ORA. However, this is not exactly the same case because even the CVE and CH units of measurement are different. Moreover, the variable, nonlinear viscoelasticity of the medium (as is the case with the cornea) can be described by different parameters, so that CVE and CH are not exactly the same. It may be that CVE is better at distinguishing one specific eye pathology and CH is better at distinguishing other pathologies. In the next step, it would be interesting to compare the CVE and CH of the same corneas, in individuals who are both healthy and with pathologies (especially glaucoma).
After examination of the characteristic values for given groups of patients, the proposed parameters can become an additional tool for ophthalmic diagnosis, which will be simple, fast, and non-invasive. In addition, the presented method of buckling determination may be useful in the creation and validation of a mechanical model of the cornea.

Conclusions
According to the best knowledge of the authors, the presented results of in vivo experiments on the eyes of healthy volunteers show for the first time that the corneal buckling effect can be observed and measured on the Corvis ST tonometer. Fine numerical analysis of the corneal profile obtained from each of the 140 video frames and examination of the dynamics of local slope variations of the profile enables a quantitative description of both buckling time moments. Some new parameters were defined and proposed regarding the processes of both bucklings and applanations. In particular, a specific parameter CVE was proposed to determine the viscous properties of the cornea due to its high reproducibility (ICC = 0.82), dimensionlessness, and high correlations with other parameters related to corneal viscoelasticity. Differences between the magnitude of this parameter for different corneas can be a useful indicator of the characteristic biomechanical feature of an individual cornea. The paper may be a prelude to the application of the proposed analysis to pathological corneas and, more generally, to an investigation of the corneal viscous properties.  Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.